A simple construction of initial data for multiple black holes 
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We consider the initial data problem for several black holes in vacuum with arbitrary momenta and 
spins on a three space with punctures. We compactify the internal asymptotically flat regions to 
obtain a computational domain without inner boundaries. When treated numerically, this leads to 
a significant simplification over the conventional approach which is based on throats and isometry 
conditions. In this new setting it is possible to obtain existence and uniqueness of solutions to the 
Hamiltonian constraint. 
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Binary black hole spacetimes are one of the great chal- 
lenges for numerical general relativity, even if no mat- 
ter sources are present. Here we consider the problem 
of finding initial data for several black holes in vacuum 
with arbitrary momenta and spins. In general relativity, 
initial data on a hypersurface cannot be specified freely, 
because the Einstein equations give rise to four equa- 
tions, three momentum constraints and the Hamiltonian 
constraint, that the initial data has to satisfy. The pur- 
pose of this note is to introduce a novel approach which is 
significantly simpler than the conventional method based 
on throats and conformal imaging. 

In all that follows we will assume vacuum, that the 
metric is conformally flat, and that the extrinsic curva- 
ture is tracefree. A convenient form of the constraints of 
general relativity can be obtained by rescaling the phys- 
ical three-metric gtu and its extrinsic curvature K? b by 
a conformal factor ip, 
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The momentum constraint becomes 
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and the Hamiltonian constraint becomes an elliptic equa- 
tion for the scalar field ip, 
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where the covariant derivatives are defined by the flat 
metric g a b, which is also used to raise and lower indices 
(see §§). 

In order to obtain black hole vacuum data, one has to 
introduce a non-trivial topology. The first calculations 
were performed by Einstein and Rosen Q in their work 
on point particles in general relativity. Various construc- 
tions for black holes based on Einstein- Rosen bridges and 
"wormholes" were given in e.g. M UJ. The spatial slice 



typically consists of two or more copies of R 3 with sev- 
eral spheres removed and identifications of the various 
spherical inner boundaries. In this way several asymp- 
totically flat regions are obtained that are connected by 
bridges or "throats" . 

The simplest example derives from the Schwarzschild 
spacetime in quasi isotropic coordinates. Considered as 
a problem on R 3 minus the point r = 0, the constraint 
equations (||) and (|J) are solved by 



m 



K ab = 0, 
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where m is the mass and r the isotropic radius. To make 
contact with the throat picture, recall that there exists 
an isometry given by r <-» m 2 /(4r) which leaves the co- 
ordinate sphere r = m/2 invariant and which maps the 
entire exterior asymptotically flat space into that sphere. 
Consequently, there exists a second asymptotically flat 
region near r = 0. Equivalcntly, one can represent this 
solution to the constraints on a space consisting of two 
copies of R 3 with a sphere excised and appropriate iden- 
tification at the spheres. 

For N black holes and non-vanishing extrinsic curva- 
ture, York and others || have developed a sophisticated 
method to solve the constraints for two asymptotically 
flat spaces that are connected by as many throats (i.e. 
excised spheres) as there are black holes, and that are iso- 
metric copies of each other. Note that there are explicit 
solutions to the momentum constraint (j^) that charac- 
terize a single black hole with given momentum P a , and 
spin S a . For example, 

3 : {P a n b + P b n a - {g ab - n a n b )P c n c ) 



K ab 
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+ ^-(e acd S c n d n b + e bcd S c n d n a ), 
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where n a is the radial normal vector. Since the conformal 
metric is flat, we will use in what follows either standard 
spherical or Cartesian coordinates, i.e. r = (x 2 + y 2 + 
z 2 ) 1 / 2 and n a — x a jr. 
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By the method of images it is possible to obtain an 
infinite series based on (0) for K ab which solves the mo- 
mentum constraint and satisfies an isometry condition at 
any number of spheres |6|^] . Given such a solution, what 
remains to be done is to solve the Hamiltonian constraint 
(||), which is an elliptic equation on R 3 minus several 
spheres, with the inner boundary given by the isometry 
and the outer boundary determined by asymptotic flat- 
ness. 

Three independent numerical implementations of the 
above scheme have been given and compared in the 
definitive paper on three dimensional black hole initial 
data by Cook, Choptuik, Dubai, Klasky, Matzner, and 
Oliveira ||. For the Hamiltonian constraint, they con- 
sider a non-linear block full approximation storage multi- 
grid scheme for Cadez coordinates, a successive over- 
relaxation scheme in Cartesian coordinates, and a mul- 
tiquadratics approach. All these approaches are greatly 
complicated by the presence of the inner spherical bound- 
aries and the isometry condition. 

Let us return to the Schwarzschild solution of the con- 
straints on a "punctured" R 3 . As noted by Misner and 
Wheeler || , and studied in detail by Brill and Lindquist 
0, the Schwarzschild solution to the constraints gener- 
alizes trivially to N black holes for time symmetry, 

N 
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oH { l , , K ab = (6) 
j= i 2 \ r - r (i)\ 

where rnt{\ characterizes the mass of the i-th black hole 
(i.e. the ADM mass is J2 m (i)) an d ^(i) 1S the location 
of the i-th black hole. For regularity of the conformal 
factor, the MWBL solution is considered on a single R 3 
with the points r — rv^ removed. We refer to the ru-\ 
as punctures. The isometry present in the Schwarzschild 
solution is lost, although there still exist minimal surfaces 
characterizing the throats . 

Let us now discuss the method that we propose to find 
data for multiple black holes with arbitrary boosts and 
spins. The idea is to compactify the internal asymptot- 
ically flat regions in order to obtain a simple domain of 
integration. Such compactification brings up issues of 
regularity of the fields (e.g. Jl(],[ll]] ) , which we address 
below. 

As before, we consider vacuum spacctimes, the metric 
is conformally flat, and the extrinsic curvature is trace- 
free. As the spatial slice we choose a single R 3 with N 
punctures as in the MWBL data. First, we solve the 
momentum constraint @ by setting 

N 

Rab = ^2 K PS(i)' (7) 

i=l 

where each term is defined by (j^) with its own origin 
r/j), momentum Pu), and spin S^y These parameters 
correspond to the ADM quantities in the limit that the 



separation of the holes is very large. The equation (0) 
defines the solution to the momentum constraint that 
we actually use, and is not just the starting point for 
the method of images that is usually invoked to distort 
K a t, to obtain an isometric solution (cmp. [fL2|, where the 
same simplification arises for a trapped surface boundary 
condition at the inner boundary). 

Given K ab as defined in (|7|), we proceed to solve the 
Hamiltonian constraint, (||). We rewrite the conformal 
factor in terms of functions a and u given by 

* — i = E^T' (8) 

a 01 l~i 2 \ r - r {i)\ 

On the punctured R , the Laplacian of 1/a is zero, so 
that the Hamiltonian constraint equation becomes (cmp. 
|llf for a single asymptotic region and vanishing linear 
momentum) 

Au + /3(1 + auy 7 = 0, (9) 
= ^a 7 K ab K ab . (10) 

To complete the definition of the problem, we have to 
specify boundary conditions for u. For asymptotic flat- 
ness at infinity we require u — 1 = 0(r~ l ) for large dis- 
tances to the punctures. 

The key question that remains is what condition we 
want to impose on u close to the punctures. As it turns 
out, to build in asymptotically flat regions as are present 
in the MWBL data near the punctures, it suffices to solve 
everywhere on R 3 without any points excised. This 
completes the statement of our proposal. 

Let us discuss existence and uniqueness of solutions to 
the modified Hamiltonian constraint equation (^) on R 3 . 
Since in this case the topology is trivial, we can show 
existence and uniqueness of a C 2 solution by repeating 
the proof given in |l| for the conventional Hamiltonian 
constraint equation on R 3 . By definition, both a and 
(3 are proportional to \r — fVj)| near the former punc- 
tures and are therefore C°, despite the fact that K a i,K ab 
goes as \f — f^\~ e at these points. From the maximum 
principle and the outer boundary condition we obtain 
that there exists at most one C 2 solution, and in par- 
ticular, that u > 1 and 1 + au > 1. Referring to |l| 
for the definition of weighted Sobolev spaces, for the ex- 
istence part of the proof we require that u G s and 
(3(1 + au)~ 7 G s+2 w ith P > 3 for the norm and 

s > 3 characterizing differentiability. The fall-off of (B) is 
the same as in the standard case. Interior points like the 
fVj) do not affect the weight S. Due to the loss of differ- 
entiability at the punctures in /3, we find that in our case 
p > 3 and s — 3. (In comparison, |l4| ] uses p — 2 but 
requires s > 4). The proof proceeds as in |l| with minor 
changes in the algebra and due care whenever it matters 
that a and j3 may vanish. So, although u is only C 2 at 
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the r*(i) , there exists a unique solution for the conformal 
factor ip = u + X/a on the punctured R 3 determined by 
our proposal for u on the unpuctured R 3 . 

Given a solution u, we can demonstrate that each 
puncture represents the "point at infinity" for another 
asymptotically flat spacetime (cmp. Hence, solving 
(^) on R 3 involves a particular compactification of N out 
of N + 1 asymptotically flat regions, one of which is dis- 
tinguished by our choice of K a b, see below. We perform a 
coordinate inversion through a sphere near the i-th punc- 
ture, f = a 2 /r, under which the metric transforms as 



V> 4 (dr 2 + r 2 dQ 2 ) = ?/> 4 {df 2 + f 2 dQ 2 ) , (11) 



with tp = tpr/a (which is not an isometry). Setting a = 
77i(j)/2 where is the bare mass of the puncture we 
are considering, we obtain for our choice of tp, (||), that 
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Therefore, the metric becomes flat as we approach the 
punctures. 

To show asymptotic flatness at the punctures, it re- 



mains to be shown that K a b = 0(f 



In the phys- 



ical variables, we just have the coordinate transforma- 
tion f = a 2 /r, = dx b /dx c = (a 2 /f 2 )L b where 
L h c = 5 b — 2n b n c . In the unphysical variables, we also 
have to take into account the transformation of ip. With 



K ab = iP 2 K p a % and K% = (AK ph ) ah , we obtain 
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Therefore, for our choice of if a 6, (Hh the momentum term 
of order r~ 2 is mapped to order f -4 while the spin term 
of order r~ 3 is mapped to order f~ 3 . This observation 
extends to multiple hole data with K a b defined in (0). 
Hence the black hole puncture data is asymptotically flat 
at the punctures. 

In fact, we also learn that the black hole seen from 
the region near the puncture appears to be unboosted, 
because the boost term is obtained from the r~ 4 term 
in the untransformed space. This term is always zero in 
our space by construction. In the conventional method, 
imposing the isometry leads to terms which go as r -4 in 
K a b. These terms may not be included in our method 
(without modification), because j3 would not be regular 
at the punctures if they were. 

Since the initial data is asymptotic to Schwarzschild 
data near the punctures, and since small spherical sur- 
faces centered at the puncture of Schwarzschild data 
are outer-trapped surfaces (with outside referring to the 
asymptotic region away from the punctures), this indi- 
cates that the region near the punctures are inside black 



holes. As in 0, several punctures may be hidden behind 
a common horizon, which we confirmed numerically but 
will not discuss here. 

We now come to the numerical implementation of our 
proposal. We found that a full approximation storage 
multigrid method built around a non-linear Gauss-Seidel 
relaxation scheme |ll| performs very well for ^ on a 
finite Cartesian grid with a Robin boundary condition. 
Note that the standard representation of Am by centered 
finite differences is only first order at the punctures. One 
could use a different prescription at the puncture points, 
but in numerical tests we found that the correspond- 
ingly lower rate of convergence is contained in a surpris- 
ingly small neighborhood of the puncture points, whether 
those are part of the grid or not. 

The multigrid program is part of BAM, a bifunctional 
adaptive mesh package for elliptic and hyperbolic prob- 
lems in three-dimensional numerical relativity (see |l6| l 
for the hyperbolic part; we do not use adaptivity here). 
Since in general the Hamiltonian constraint has to be 
solved numerically, it is important to provide at least 
one efficient numerical implementation. There is no rea- 
son to believe that the multigrid method is the only good 
method. However, when available, multigrid methods 
have proven to be among the best performers. We con- 
sider the absence of irregular boundaries in our method 
to be a valuable feature since it makes a straightforward 
multigrid implementation possible. 

As a test, let us compute the correction to the confor- 
mal factor for two equal mass black holes boosted toward 
one another when we keep only terms of small P (the 
momentum of each hole) and L (the coordinate positions 
of the holes along the z-axis being ±L). Consider (3 to 
be of order e, and consider the MWBL solution to be 
the zeroth order solution to the Hamiltonian constraint, 



t(0) 



t(i) 



1 + We now have 



eAu(i) = -/3(l + a) 
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The ADM mass to first order in e is 
M ADM = -^-i Va(-) dS a + 1- 



P{l + a)- r dV 
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This equation predicts the correction to the ADM mass 
resulting from the full non-linear solve to within about 
5% when P/m < 1.0 and L/m < 0.2 and the grid is a 
cube whose sides have a length of about 16m and are 
resolved with 67 zones. 

This calculation provides one test for the correctness 
of the solution of the data set, at least for the ADM 
mass of the spacetime for perturbative cases. This is 
valuable, because the code cannot be tested for the case 
when K a b = since then the solution is trivially the 
MWBL solution. 
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A2B8, box of size 32 
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4.000 


ii 3 


1.159e-08 


3.813e-03 


4.029 


9.209 


2.000 


19 3 


5.350e-08 


1.984e-03 


-2.119 


10.097 


1.000 


35 3 


2.141e-07 


5.400e-03 


0.838 


10.601 


0.500 


67 3 


2.582e-07 


1.216e-02 


1.121 


10.833 


0.250 


131 3 


1.552e-07 


2.222e-02 


1.430 


10.917 


0.125 


259 3 


6.629e-08 


3.383e-02 


1.712 


10.942 



TABLE I. Results for the A2B8 data set of [9]. All units 
are in terms of m( 2 ), h is the grid spacing, n the number of 
grid points, d is the residual, r the truncation error, a the 
convergence rate, and Madm is the ADM mass of the box. 



As a general 3d strong field example we present in Tab. 
| data comparable to the A2B8 data set of ||. Note 
that the data sets produced by the two methods are not 
identically the same. They are equally generic, but each 
has a different gravitational wave content. 

There are two punctures with different mass parame- 
ters of order one, and with comparatively large values for 
boosts and spins. The results are obtained for a single 
box centered around the origin for various grid spacings. 
Each run was performed for the same multigrid param- 
eters such that the values for the residua reflect the ef- 
ficiency with which the discretized equations are solved 
at different resolutions. We chose to make the residua 
much smaller than the truncation error estimate, which 
turns out to be in a range compatible with the still rather 
large grid spacing. The convergence rates a are based on 
u at three resolutions. As expected, we found in this and 
other examples that the convergence rate at the punc- 
ture points was one order less than elsewhere, but we 
also found that convergence was only affected very close 
to the punctures. As the resolution increases, the mass 
estimate converges with a rate comparable to a. 

To summarize, we have introduced a new setup for the 
initial data problem for multiple black holes with arbi- 
trary mass, momentum, and spin in three dimensions. 
Comparing with the wormhole constructions of black 
holes, the key difference is how the various asymptotically 
fiat regions are defined on a single copy of R 3 . Instead 
of using an isometry condition at interior spheres, we re- 
move the inner boundary by compactification and look 
for a solution to a modified Hamiltonian constraint equa- 
tion on R 3 , for which we have existence and uniqueness. 
A corresponding proof should be possible for the isomet- 
ric wormhole data, although to our knowledge, such a 
proof as not been given yet. 

There are two obvious numerical advantages compared 
to the wormhole approach. First note that the solution 
K ab to the momentum constraint is all we need, i.e. since 



we do not impose the isometry condition there is no need 
for the method of images. Second, in the solution of 
the Hamiltonian constraint we avoid the numerical com- 
plications due to an inner boundary. Compared to the 
standard Cartesian (Cadez) method, one (two) levels of 
sophistication less are required for an implementation, 
which should make general black hole initial data more 
widely available. 

Finally, note that black hole puncture data has already 
been successfully evolved in the case of 3d Schwarzschild 
initial data (jl^,[l6| , which is possible even when the punc- 
tures are part of the numerical grid since there are no 
physical singularities on the initial slice. For a long term 
evolution of black holes, one may start with data ob- 
tained by either the throat or puncture construction, and 
then cut out the interior regions at the apparent horizon 
to avoid the physical singularities on future slices Jl2],[l8| . 
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